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An implicit, fully characteristic, numerical scheme for solving the field equations of a cosmic string 
coupled to gravity is described. The inclusion of null infinity as part of the numerical grid allows 
us to apply suitable boundary conditions on the metric and matter fields to suppress unphysical 
divergent solutions. The code is tested by comparing the results with exact solutions, checking that 
static cosmic string initial data remain constant when evolved and undertaking a time dependent 
convergence analysis of the code. It is shown that the code is accurate, stable and exhibits clear 
second order convergence. The code is used to analyse the interaction between a Weber-Wheeler 
pulse of gravitational radiation with the string. The interaction causes the string to oscillate at 
frequencies proportional to the masses of the scalar and vector fields of the string. After the pulse 
has largely radiated away the string continues to ring but the oscillations slowly decay and eventually 
the variables return to their static values. 
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This is the second paper in a series devoted to the study of dynamic cosmic strings. In the previous paper [0, 
henceforth referred to as paper I, we derived the equations of motion for a time dependent cylindrically symmetric 
cosmic string coupled to a gravitational field with two degrees of freedom. The treatment involved using a Geroch 
i decomposition to reformulate the problem in terms of matter fields on a 2 + 1 dimensional spacetime and two geo- 
metrical variables v and r which describe the gravitational degrees of freedom. Unlike the original cylindrical metric 
the reduced 2 + 1 spacetime is asymptotically flat and this allows us to implement a numerical treatment of the field 
equations which includes null infinity as part of the grid. In paper I we presented a Cauchy-characteristic matching 
(CCM) code that reproduced the results of several vacuum solutions with both one and two degrees of freedom but 
O" 1 did not perform entirely satisfactorily in the presence of matter. We have therefore developed a second implicit, 
| \ purely characteristic code. The numerical details, the testing and the convergence analysis of this implicit code will 
bJ[). be described in this paper. We also give a detailed analysis of the interaction between the string and a pulse of 
' gravitational radiation. 

• i-H . After briefly describing the variables used for the coupled system in section II the details of the numerical scheme 
are presented in section III. We begin by describing the relaxation method we use for the much simpler problem of a 
static cosmic string in Minkowski spacetime and demonstrate how this approach leads naturally to the implicit scheme 

- - ■ used to solve the general case of a dynamical cosmic string coupled to gravity. The testing of the code is described 
in section IV. This involves comparing it with exact solutions, checking that static cosmic string initial data remain 
static when they are evolved, and undertaking a time dependent convergence analysis of the code. We are able to show 
that the code is accurate, stable and shows clear second order convergence. In section V we analyse the interaction 
between an initially static cosmic string and a Weber-Wheeler type pulse of gravitational radiation. We show that 
this interaction causes the string fields S and P to oscillate and examine how this oscillation depends upon both the 
width and amplitude of the pulse and also upon the constants e, A and rj. The key result is that the frequencies are 
essentially independent of the nature of the Weber- Wheeler wave but are proportional to the masses rag and mp 
associated with the scalar field S and the vector field P. The string continues to ring after the gravitational pulse has 
largely radiated away, but the oscillations slowly decay and the variables eventually return to their static values. It 
is only possible to observe this effect because of the long term stability of the code. Finally in section VI we discuss 
our results and outline future work. 
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II. THE VARIABLES AND FIELD EQUATIONS 



We begin with a brief summary of the variables used to describe the spacetime and the string. In cylindrical polar 
coordinates (i, p, </>, z) we may write the line element for a cylindrically symmetric spacetime with two degrees of 
freedom as a modified version of that given by Jordan, Ehlers, Kundt and Kompaneets 0, p] 



dg 2 = e 2(7"V)(^2 _ dp2) _ p 2 e -24, d4> 2 _ e 2W+/0 ( w< ty + dz f 



(1) 



where tp, to, p and 7 are functions of i, p. As shown in paper I, in order to work with an asymptotically flat spacetime 
it is convenient to make a Geroch decomposition W in which the 4-dimensional metric is replaced by a 3-dimensional 
one and two auxiliary scalar fields. These scalar fields are the norm of the axial Killing vector v and the Geroch 
potential r and are related to ip and to by equations (15) and (16) of paper I, while the conformal 3-dimensional line 
element is given by 



da 2 = e 2 ^(dt 2 - dp 2 ) - p 2 e 2 »dcj> 2 . 



(2) 



The cosmic string is described by a complex scalar field $ coupled to a U(l) gauge field A^, which in cylindrical 
symmetry can be written as fy] 



A lt = -[P(t,p)-l]V li <i 



(3) 
(4) 



It is also helpful to introduce rescaled coupling constants and variables given by 

S 



X 



(5) 



„- - , 

r = V~\r)p, 
t = VXrjt. 



(6) 

(7) 
(8) 



Here rj is the vacuum expectation value of the scalar field while a represents the relative strength of the coupling 
between the scalar and vector field given by e, compared to the self-coupling of the scalar field given by A. Critical 
coupling, for which the masses of the scalar and vector fields are equal, is given by a — 8 || . 

The field equations for a cosmic string coupled to gravity were derived in paper I in terms of Cauchy coordinates 
[Paper I (26)~(33)] and also in terms of compactified characteristic coordinates [Paper I (43)-(49)]. However the fully 
characteristic numerical scheme also requires characteristic equations in the inner region. In terms of the retarded 
time u = t — r and the radius r the field equations are given by 
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where □ represents the flat-space d'Alembert operator 

_ d 2 d 2 1 / d d 
dudr dr 2 r \dr du 

Note that there are two further Einstein equations which are not used in the numerical scheme as they are a conse- 
quence of the above equations and their derivatives. These equations are only used to provide a check on the numerical 
accuracy of the code. 




III. NUMERICAL METHODS 



In order to solve the above field equations we have developed two independent codes. The first is based on the 
Cauchy-characteristic matching code of Dubai et al. Q| and d'Inverno et al. g]. This code performs well in the 
absence of matter and has been used in paper I to study several cylindrically symmetric vacuum solutions. In paper 
I we have also given details of the convergence analysis and described the modifications with respect to B that lead 
to long term stability with both gravitational degrees of freedom present. However the CCM code performed less 
satisfactorily in the evolution of the cosmic string. This is due to the existence of unphysical solutions to the evolution 
equations (|9|)-(|l4]) which diverge exponentially as r — ► oo. Controlling the time evolution near null infinity by means 
of a bump function enabled us to select the physical solutions with regular behavior at / + , but the bump function 
itself introduced noise which eventually gave rise to instabilities. We therefore implemented a second implicit, purely 
characteristic, code which allowed us to directly apply boundary conditions at the origin as well as null infinity and 
thus suppress diverging solutions. It is interesting that this problem is already present in the calculation of the static 
cosmic string in Minkowski spacetime. We will, therefore, first describe the numerical scheme used in the static 
Minkowskian case where the equations are fairly simple. We then present the modifications necessary for the static 
and dynamic case coupled to the gravitational field. 



A. The static cosmic string in Minkowski spacetime 

In we set the metric variables to their Minkowskian values and all time derivatives to zero to obtain the 

equations for the static cosmic string in Minkowski spacetime (cf. |^|) 

r-ifr- 1 ^) =«X*P, (16) 



dr \ dr 



d ( dX\ 
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dr V dr J 



X [P 2 + Ar 2 {X 2 - 1)] . (17) 



The boundary conditions are 



P(0) = 1, lim P(r) = 0, 

r—*OG 

X(0) = 0, lim X{r) = 1. (18) 

r — >oc 

In order to cover the whole spacetime with a finite coordinate range, we divide the computational domain into 
two regions. In the inner region (0 < r < 1) we use the coordinate r, while in the outer region we introduce the 
compactified radius 

y = ± (19) 

which covers the range 1 > y > corresponding to the region r > 1. It is also useful to combine r and y into a single 
radial variable w defined by 

r for < r < 1 , . 

3 - for r > 1, {ZU > 
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FIG. 1. The combined grid of the inner and the outer region. Note that both grid points, n\ and ni + 1, correspond to the 
position r = 1 <4> y — 1. These points form the interface of the code and facilitate transformation of the variables from the 
coordinate system r into that of y. 



so that the region < w < 3 corresponds to < r < oo and infinity is mapped to w = 3 Q . 
In terms of the variable y, and (^) take the form 
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(21) 
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The number of grid points in each region may differ, but each half-grid is uniform. Thus we use a total of N := n\ +112 
grid points where the points labelled n\ and n\ + 1 both correspond to the position r = 1 = y. The points ni, ni + 1 
form the interface between the two regions (see Figure ^). One point will contain the variables in terms of r, the 
other in terms of y. With the computational grid covering the whole spacetime, we now face a two point boundary 
value problem. Due to the existence of unphysical solutions diverging at y = shooting methods turned out to be 
unsuitable for solving this problem. On the other hand numerical relaxation, as described in |J for example, allows 
us to directly control the behavior of P and X at infinity. The form of equations ([l6]), ( |l7| ) suggests that in order 

fn iTtrT-no a fivcf r\rr\or ctrc+om -roro cVinnlrl inf rr\r\ nco -flio cinvilicinr -\ra ri a VJ oc ( ) — v P, aud R — vX The 



to write them as a first order system we should introduce the auxiliary variables Q — r 
equations may then be written in the form 

P.r = rQ, 



(23) 
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The corresponding equations in the outer region are given by 
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Standard second order centered finite differencing results in 4 (AT — 2) non-linear algebraic equations which are sup- 
plemented by the 4 boundary conditions (18) and 4 interface relations 



(31) 
(32) 
(33) 
(34) 

We then start with piecewise linear initial guesses for P and X (and the corresponding derivatives Q and R) and 
solve the AN algebraic equations iteratively with a Newton-Raphson method. Results for various choices of the string 
parameter a are shown in Figure 6 of paper I. 
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B. The static cosmic string coupled to gravity 

From the numerical point of view, the problem of solving for a static cosmic string coupled to gravity through the 
Einstein equations is virtually identical to that of a static string in Minkowski spacctime. The only difference is the 
much higher degree of complexity of the equations due to the appearance of the functions v, r, /i and 7 as extra 
variables. We do not present the equations here, since they may be derived from the fully dynamic case by setting 
all time derivatives to zero in the relevant equations. The solution is again obtained using the relaxation method 
described in the previous section. As our initial guess for the metric variables we use Minkowskian values, and for 
the string variables X and P we use the previously calculated values for a Minkowskian string with the same string 
parameters. The results are shown in Figure 7 of paper I. 



C. The dynamic code 

In the dynamic case all variables v, r, /1, 7, P and X are functions of u, r and we have to solve the system (|^)- 
( p"4| ) of partial differential equations (PDEs). In order to control the behavior of the solution at infinity, we need 
a generalisation for PDEs of the relaxation scheme applied to ordinary differential equations (ODEs). In view of 
the characteristic feature of the relaxation scheme, namely the simultaneous calculation of new function values at all 
grid points, this generalisation leads directly to implicit evolution schemes as used for hyperbolic or parabolic PDEs. 
Therefore, the dynamic code is based on the implicit, second order in space and time Crank-Nicholson scheme (see ^ 
for example). For each solution step we consider two spatial slices of the grid- type of Figure [IJ labelled n and n + 1. 
We then apply centered finite differencing according to 

fn+1 1 fn 

f= h+1 +h ' Jk+1 ' Jfc , (35) 

fn+1 fn+1 1 rn _ fn 

f,r = U+1 \ A h+1 h , (36) 



where / represents any of our variables. Assuming that all functions are known on slice n we arrive at a large set of 
algebraic equations for the similar to the static case, which needs to be supplemented by boundary conditions 

and interface relations analogous to ([5l])-(|34|). Again the system of algebraic equations is solved iteratively with the 
Newton-Raphson method. The initial guess for the data on the new slice n + 1 is taken from the previous slice and 
convergence is typically achieved within three iterations. For this purpose we rewrite the dynamic equations (p[)-(p^) 
as a first order system. The equations for the variables v, r and X involve radial derivatives which may be written 
in terms of the second order operator ■§^{t-§^) so we introduce the corresponding variables A^ = rz^ r , T = rr r and 
R = rX_ r [cf. paper I equations (110)— (115)] . The equation for /x on the other hand involves ■§p( 7 ' 2 -§^), while that for 
P involves ■§^( r ~ 1 -§p)- We therefore introduce the corresponding variables M — r 2 /i r and Q = r~ 1 P_ r . Finally the 
equation for 7 involves only one r derivative and may be written in first order form without having to introduce any 
further quantities. In terms of these variables (p|)-([T4|) become: 
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The corresponding first order system in the outer region is given by 
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In order to solve these equations we must supplement them by appropriate initial and boundary conditions. We start 
by considering boundary conditions on the axis. In general we find the code is more stable if one imposes boundary 
conditions on the radial derivatives rather than the variables themselves. For the variables v, t and X we therefore 
impose the required boundary conditions on the initial data, but in the subsequent evolution we impose the weaker 
condition that their radial derivatives are finite on the axis. This ensures that the evolution equations propagate the 
axial conditions given on the initial data. For the variable \x we impose the condition that M is zero on the axis 
which is equivalent to the rather weak condition that r 2 /u r vanishes there. The inverse power of r in the definition of 
Q makes it unsuitable to specify the value of this quantity at r — so in this case we work with the variable directly 
and require that P — 1 on the axis. Finally the variable 7 is given by a purely radial equation, so in this case we must 
specify the value on the axis which is chosen to be zero to ensure elementary flatness. Therefore at r — we require 



N = 0, 
T = 0, 
M = 0, 
7 = 0, 

P = l, 
R = 0. 



(60) 
(61) 
(62) 
(63) 
(64) 
(65) 



For the boundary conditions at null infinity we know that regular solutions of the cylindrical wave equation have 
radial derivatives that decay faster than 1/r so that we may take the variables N, T and R, which satisfy a wave type 
equation, to vanish at y = 0. The asymptotics of [i are slightly different due to the additional power of r in the radial 
derivative (similar to the spherically symmetric wave equation) but for a regular solution fj, tV vanishes at null infinity. 
The P equation does not satisfy a wave type equation due to the inverse power of r but has asymptotic behavior given 
by a modified Bessel function. The physically relevant finite solution has exponential decay so in this case one may 
impose the condition that Q = at y — 0. Hence we require the solution to satisfy the following boundary conditions 
at y = 



AT = 0, 
T = 0, 

= °> 
Q = 0, 
R = 0. 



(66) 
(67) 
(68) 
(69) 
(70) 



These boundary conditions are sufficient to determine the solution of the first order system (|38D-(p9|) while suppressing 
the unphysical solutions which are singular on the axis or null infinity. Note that 7 is determined by the constraint 
equation (JL2h, which is a first order ODE, and thus only needs one boundary condition. 
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IV. TESTING THE CODE 

In paper I we used the CCM code to reproduce several exact vacuum solutions to an accuracy of about 10~ 5 with 
second order convergence. We have also checked the properties of the codes for the static cosmic string in Minkowski 
spacetime and coupled to gravity. Both clearly showed second order convergence. Here we will focus on testing the 
implicit dynamic code. We have carried out four different independent tests, namely 

(1) Reproducing the non-rotating vacuum solution of Weber and Wheeler p0| , 

(2) Reproducing the rotating vacuum solution of Xanthopoulos @ , 

(3) Using the results for the static cosmic string (paper I) as initial data and checking that the system 
stays in its static configuration, 

(4) Convergence analysis for the string hit by a Weber-Wheeler wave. 

Two additional tests arise in a natural way from the field equations and the numerical scheme. As described above 
there are two additional field equations which are algebraic consequences of the other field equations. We have verified 
that these equations are satisfied to second order accuracy (~ Ar 2 ). Furthermore the numerical scheme calculates the 
residuals of the algebraic equations to be solved, which have thus been monitored in test runs. They are satisfied to a 
much higher accuracy (double precision machine accuracy), so the total error is dominated by the truncation error of 
the second order differencing scheme. Another independent test is the comparison with the explicit CCM code which 
yields good agreement for as long as the latter remains stable. The four main tests are now described in more detail. 



A. The Weber— Wheeler wave 



In the first test we evolve the analytic solution given by Weber and Wheeler |1C(], which describes a gravitational 
pulse of the + polarisation mode. This solution has two free parameters, a and 6, which can be interpreted as the 
width and amplitude of the pulse. The equations together with a more detailed discussion have been given in paper I. 
We prescribe v as initial data according to the analytic expressions obtained for a = 2 and b = 0.5 and set the other- 
free variables to zero, while 7 is calculated via quadrature from the constraint equation (|l^). In Figure || we show 
the deviation of the numerical results from the analytical one for N = 1920 grid points (320 points in the inner, 1600 
points in the outer region) and a Courant factor of 0.5 with respect to the inner region. The convergence analysis 
(see below) shows that this number of points provides sufficient resolution in the outer region while still keeping 
computation times at a tolerable level. All computations presented in this work have been obtained with these grid 
parameters, unless stated otherwise. The code stays stable for much longer time intervals than shown in the figure, 
but does not reveal any further interesting features as the analytic solution approaches its Minkowskian values and 
the error goes to zero. 
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B. The rotating solution of Xanthopoulos 



Xanthopoulos JTl] ] derived an analytic vacuum solution for Einstein's field equations in cylindrical symmetry con- 
taining both the + and x polarisation mode. Its analytic form in terms of our metric variables and a more detailed 
discussion has been given in paper I. The solution has one free parameter a which is set to one in this calculation. 
The error of our numerical results is displayed in Figure [?], where we have used the same grid parameters as in the 
Weber-Wheeler case. Again we have run the code for longer times and found that the error approaches zero. We 
conclude that the code reproduces both analytic vacuum solutions with excellent accuracy comparable to that of the 
CCM code and exhibits long term stability. 



C. Evolution of the static cosmic string 



The tests described above only involve vacuum solutions, so the matter part of the code and the interaction between 
matter and geometry has not been tested. An obvious test involving matter and geometry is to use the result for 
the static cosmic string in curved spacetime as initial data and evolve this scenario. All variables should, of course, 
remain at their initial values. We have evolved the static string data for our standard grid and the parameter set, 
a = 1 and rj = 0.2, which corresponds to a strong back-reaction of the string on the metric. The results are shown in 
Figure 0. The system stays in its static configuration with high accuracy over a long time interval. 



D. Convergence analysis 

Our investigation of the interaction between the cosmic string and gravitational waves will focus on the string being 
hit by a wave of the Weber- Wheeler type. In order to check this scenario for convergence we have run the code for 
the parameter set r\ — 0.2, a = 1, a = 2, b = 0.5 for different grid resolutions, where a and b are again the width 
and amplitude of the Weber- Wheeler wave. In our case it is of particular interest to investigate the time dependence 
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FIG. 4. The deviation of the metric and matter variables from the initial data in the case of evolving a static cosmic string 
with a = 1, T] = 0.2. For our standard grid with 1920 points, the configuration stays static to an accuracy of about 10 over 
a range of more than 30000 time steps. 

of the convergence to see what resolution we need in order to obtain reliable results for long runs. We calculate the 
convergence rate in the same way as in the static case (cf. paper I), but this time the £ 2 _n orm is a function of time. 
So for each variable we have 

, 3[A *» )W . y SgWEM! , (71) 

where the upper label "4320" indicates that the high resolution reference solution has been calculated for N = 4320 
grid points. In Figure || we show the convergence factor ^I 1 ^ 1920 ]/^! 1 ^ 2880 ] as a function of u for v, fi, 7, P and 
X. The initial data for r is identically zero for this scenario and stays zero during the evolution. The number of 
grid points is increased by a factor of 1.5 here (instead of the more commonly used 2) to reduce the computation 
time. Only points common to all grids have been used in the sum in equation (fzT|) . For second order convergence we 
would expect a convergence factor of 1.5 2 . Although the results in Figure [5] show weak variations with u, second order 
convergence is clearly maintained for long runs. In each case the outer region contains 5 times as many grid points 
as the inner region (e.g. n\ — 320, rc.2 = 1600 for the N — 1920 case). The reason for this is that in the dynamical 
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FIG. 5. The convergence factor £ 2 [^/ 1920 ]/£2[ } if 2SS0 ] is plotted as a function of u. We expect a convergence factor of 2.25 
since the number of grid points is multiplied by 1.5. Even though our results show weak variability at later times, second order 
convergence is maintained throughout long runs (more than 30000 time steps with N = 1920). 

evolutions X and especially P exhibit significant spatial variations out to large radii. Due to the compactification, 
the spatial resolution of our grid decreases as we move towards null infinity and to resolve the spatial variations of the 
string variables out to sufficiently large radii we therefore have to introduce a large number of grid points in the outer 
region. No such problems occur in the inner region. If significantly fewer grid points are used in the outer region for 
this analysis, the convergence properties of the string variables can deteriorate to roughly first order level. 



V. TIME DEPENDENCE OF THE STRING VARIABLES 
A. Static cosmic strings excited by gravitational waves 



The scenario we are going to investigate in this section is an initially static cosmic string hit by a gravitational wave 
of Weber- Wheeler type. For this purpose we use the static results with two modifications as initial data. Firstly the 
static metric function vq is multiplied by the exact Weber- Wheeler solution to simulate the gravitational wave pulse. 
Thus we guarantee that initially the cosmic string is indeed in an equilibrium configuration provided the wave pulse 
is located sufficiently far away from the origin and its interaction with the string is negligible. Ideally the numerical 
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FIG. 6. The initial data for v, P and X at uo = —20 for the standard parameters a = 1, r} = 10 -3 , a — 2, b = 0.5. The 
gravitational wave pulse is located in a region where the string fields P and X have almost fallen off to their asymptotic values. 



calculation would start with the incoming wave located at past null infinity. In order to approximate this scenario, we 
found it was sufficient to use the large negative initial time uq = —20. The second modification is to calculate 7 from 

) to preserve consistency with the Einstein field equations. In Figure ^ the corresponding 
initial data for v, P and X are shown for the parameter set 77 = 10~ 3 , a = 1, a = 2 and b = 0.5. From now on we 
will refer to these values as "standard parameters" and only specify parameters if they take on non-standard values. 
Note that r vanishes on the initial slice in this case and stays identically zero throughout the evolution. The case of 
rotating gravitational waves hitting a cosmic string will be analysed in a future publication. The time evolution of 
the "standard configuration" is shown in Figure ^ where we plot v 7 /x, 7, P and X as functions of w at times —20, 
— 10, 0, 2 and 10. While the wave pulse is still far away from the origin, its interaction with the cosmic string is 
negligible (dotted lines). When it reaches the core region, however, it excites the cosmic string and the scalar and 
vector field start oscillating (dashed curves). After being reflected at the origin, the wave pulse travels along the 
outgoing characteristics and the metric variables v, fi and 7 quickly settle down into their static configuration which 
is close to Minkowskian values for rj — 10~ 3 . The vector and scalar field of the cosmic string, on the other hand, 
continue ringing albeit with a different character. Whereas the oscillations of the scalar field X are dominant in the 
range r < 2 and have significantly decayed at u — 10 as shown in the figure, the vector field oscillations propagate to 
large radii and fall off very slowly (solid curves). This behavior is also illustrated in the right panel of Figure |8| which 
shows a contour plot of P as a function of (u, r) out to r = 50. We shall see below, that the oscillations of P will also 
decay eventually and the cosmic string will asymptotically settle back into its equilibrium configuration. 



the constraint equation (12 



B. Frequency analysis 

We will now quantitatively analyse the oscillations of the scalar and vector field of the cosmic string. Since we are 
working in rescaled coordinates, time and distance are measured in units of l/VXr] and frequency in its inverse. To 
avoid complicated notation, however, we will omit the units from now on unless the meaning is unclear. In order to 
measure frequencies, we Fourier analyse the time evolution of the scalar and vector field for fixed radius r. Figure ^ 
shows P and X for standard parameters as functions of time at r = 1 together with the corresponding power spectra. 
The plots for X show a characteristic frequency fx — 0.43 whereas for P we find a strong peak at fp — 0.16. The 
spectrum for P, however, also shows a strong mode at / = 0.43. We have calculated similar power spectra for a 
large class of parameter sets and discovered this effect on numerous occasions — in addition to a strong peak at the 
characteristic frequency of P or X there is a second maximum at the characteristic frequency of the other field. In 
general the characteristic mode of X resulted in stronger peaks at smaller radius, that of P was stronger at larger 
radii. We attribute this feature to the interaction between the scalar and vector component of the cosmic string. The 
variation of the relative strength of the oscillations with radius confirms the corresponding observation in Figure [?]. 
The accuracy of the measurements of the frequencies is limited by the resolution of the Fourier spectra which again 
is limited by the time interval covered in the evolution and, thus, by computation time. For tolerable computation 
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FIG. 7. The metric and string variables are plotted as functions of w for it = —20 (dotted), u = —10 (long dotted), u = 
(dashed), u = 2 (long dashed) and u = 10 (solid line). For clarity the the graphs of P are distributed over two panels. The 
wave pulse (in v) initially moves inwards. It excites the string, is reflected at the origin and travels outwards. After u = 10 
only P differs significantly from the static configuration as the oscillations slowly decay and propagate towards larger radii (cf. 
Figure ||) . 



times we get an accuracy A/ « 0.01 which corresponds approximately to one bin in the frequency spectra. 

In order to investigate the dependency of the oscillations on a, 77, a, b and the radial position r, we have varied 
each parameter over at least two orders of magnitude while keeping the other parameters at standard values. We 
have found the following dependencies: 

(1) The frequencies of both X and P did not show any variations with rj for rj < 0.1. (Note that r\ does, 
however, appear in the units). For larger values of 77, the non- linear interaction between string and 
geometry becomes dominant and we did not find a simple relation between frequency maxima and 
parameters. 

(2) The variation of the parameters a and 6, the width and amplitude of the Weber- Wheeler pulse, 
has no measurable effect on the frequencies of X and P, but only determined the amplitude of the 
oscillations. A narrow, strong pulse leads to larger amplitudes. 
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FIG. 8. The cosmic string variable P is shown as a function of radius and time for a = 0.2 (left) and a = 1 (right) (all other 
parameters have standard values). Note that we use the radial variable r out to r — 50 here. The ringing can clearly be seen 
and shows a lower frequency for smaller a. 
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FIG. 9. Upper panels: The variables P and X at r = 1 are plotted as functions of u for a = 1, 77 = 10~ 3 , a — 2 and b = 0.5. 
Lower panels: The corresponding power spectra. The spectrum for P shows a strong peak at / = 0.16, that for X at / = 0.43. 
The latter component, however, is also clearly present in the spectrum of P which we attribute to the interaction between the 
fields (see text for details). Note that due to our rescaling of the coordinates, u is measured in units of I/VX77. 

(3) For small r the oscillations in X are stronger, whereas those for P dominate at large r. The frequency 
values, however, do not depend on the radius. For radii greater than 10 the oscillations in X had 
decayed so strongly that we could no longer measure its frequency. 

(4) To first order approximation the frequency of the scalar field, fx , is independent of a over the observed 
range 0.1 < a < 10. fp, however, shows a strong dependence on a which is also illustrated in Figure 
pi where we compare contour plots of P obtained for a = 0.2 and 1. The frequency is significantly 
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FIG. 10. The frequencies fp ('+') and fx ('x') are plotted as functions of a. The curves have been obtained from fitting 
power laws. For 5 < a < 8 the expected frequencies for P and X are similar and we found only one maximum in the power 
spectra. These values are plotted as filled lozenges and have not been included in the regression analysis as their classification 
is not obvious. The power law for P coincides with the observed values to within measurement accuracy. The frequencies of 
X, however, show a significant, albeit small, deviation from the fitted constant line. 



larger for a = 1. 

In Figure [l^ we show fp and fx as functions of a together with the power laws obtained from a linear regression of 
the corresponding double logarithmic data. In the range 5 < a < 8 (where a = 8 corresponds to the critical coupling), 
the expected values of fp and fx become similar and we observed only one maximum in the corresponding power 
spectra. Therefore a classification with respect to the vector or scalar origin of the frequencies is not obvious. These 
values (shown by filled lozenges in Figure |l0| ) have not been used in the regression analysis. We obtain power law 
indices ox — 0.00 and ap = 0.50, so that 

fx ~ const. (72) 
fp ~ sfa. (73) 

Whereas the fitted curve for fp coincides with the observed values to high accuracy, Figure |l^ shows a small but 
significant deviation of the measured fx from the fitted constant line. Indeed the sizable range over a for which 
we observe only one frequency indicates significant non-linear interaction similar to the effect of phase locking in 
non-linear systems of ODEs [jl2j. If we consider fx and fp to be given by ( |72| ) and ( |73|) in terms of the rescaled 
unphysical coordinates and we transform this back into physical units using a = e 2 /A, we arrive at the following 
relations for the physical variables 

fx ~ VA?7, (74) 



fp ~ e V . (75) 

As shown in |^] up to constant factors y/\r) and erj are the masses of the scalar and the vector field, rax and rap, and 
we conclude that X and P have characteristic frequencies 

fx ~ mx, (76) 



f P - rap. (77) 

Since the frequencies for X and P seem only to depend upon the respective masses we have attempted to confirm 
these results by considering the oscillations of a cosmic string in two further scenarios. Firstly since the frequencies 
do not depend upon the Weber- Wheeler pulse we take as initial data the static values for the metric variables but 
excite the string by adding a Gaussian perturbation to either the X or P static initial values. The evolution is 
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FIG. 11. The upper four plots show the difference between the evolved functions u, /i, 7 and X and their corresponding 
static results. For P a similar 3-dimensional plot is not suitable since it fails to resolve the oscillations of the vector field. 
Therefore we plot the ^2-norm (dashed line) and the maximum (solid line) of AP as a function of time, v, fi, 7 and X quickly 
settle down in their equilibrium configuration to numerical accuracy. The decay of the oscillations of P takes significantly more 
time but eventually P also approaches its equilibrium state. 



then computed using the fully coupled system. Secondly since the frequencies do not seem to depend upon the 
strength of the coupling to the gravitational field we have completely decoupled the gravitational field and considered 
the evolution of a cosmic string in Minkowski spacetime. The initial data is taken to be that for a static string in 
Minkowski spacetime with a Gaussian perturbation to either the X or P values. The evolution is then computed 
using the equations for a dynamical string in a Minkowskian background [equations (85) and (86) of paper I]. In both 
cases we find the same frequencies, to within an amount A/ = 0.01, that we found in the original case of the fully 
coupled system excited by a Weber-Wheeler pulse. Furthermore the frequencies did not depend on the location or 
shape of the field perturbation nor upon the choice of X or P as the perturbed field. 
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C. The long term behavior of the dynamic string 



The time evolution shown in Figures [7] and |8| indicate that the oscillations of the cosmic string excited by gravita- 
tional waves gradually decay and metric and string settle down into an equilibrium state. We have calculated a very 
long run (—20 < u < 410) to investigate the long term behavior in detail. The unphysically large value of r\ = 0.1 
is chosen for this calculation in order to guarantee a strong interaction between spacetime geometry and the cosmic 
string. In Figure [ll] we show the difference A/ := f cvo \ — / s tat between the time-dependent v, /i, 7 and X and their 
corresponding static results obtained for the same parameters. For the vector field P a similar 3-dimensional plot 
would require an extreme resolution to properly display the oscillations of the vector field (cf. Figure 0). For this 
reason we proceed differently and calculate the ^2-norm and the maximum of AP for each slice u = const. Both 
functions are plotted versus time in Figure |H} The incoming wave pulse can clearly be seen as a strong deviation 
of v from the static function. The pulse excites the cosmic string and is reflected at the origin at u = 0. The 
metric variables and the scalar field X then quickly reach their equilibrium values. The oscillations in P decay on a 
significantly longer time scale which is also evident in Figures [7] and || and the ^2-norm of AP slowly approaches zero. 
Significantly longer runs than shown here are prohibited by the required computational time, but the results indicate 
that P will also eventually reach its equilibrium configuration. 

VI. CONCLUSION 

In this paper we have described the details of the implicit, fully characteristic, numerical scheme which is used to 
solve the field equations for a cosmic string coupled to gravity. A feature of the cosmic string equations is that they 
admit exponentially diverging unphysical solutions. By using a Geroch decomposition it is possible to reformulate 
the problem in terms of fields which describe the string on an asymptotically flat 2 + 1-dimensional metric as well as 
two auxiliary fields v and r which describe the gravitational degrees of freedom. We can then introduce a conformally 
compactified radial coordinate y which allows us to include null infinity as part of the numerical grid. As well as 
avoiding the need to introduce artificial outgoing boundary conditions at the edge of the grid this approach has the 
advantage that we can also enforce boundary conditions for the string variables at null infinity which rule out the 
unphysical solutions. The use of the geometrically defined Geroch variables also improves the long term stability of 
the code compared to the use of metric variables. 

The code has been shown to reproduce the results of two exact vacuum solutions, the Weber- Wheeler solution which 
describes a pulse of gravitational radiation with just the + polarisation state, and a solution due to Xanthopoulos 
which describes a gravitational wave with both the + and x polarisation states. The code has also been shown to 
reproduce the results of the static cosmic string code in that initial data corresponding to a static solution do not 
change when evolved in time using the dynamical code. For both the exact vacuum solutions and the static initial 
data the code shows excellent long term stability. Finally a time dependent convergence analysis demonstrates clear 
second order convergence of the code. 

After demonstrating the reliability of the code we use it to analyse the interaction between an initially static 
cosmic string and a Weber-Wheeler type pulse of gravitational radiation. We find that the gravitational wave excites 
the string and causes the string variables X and P to oscillate before the configuration slowly settles back into its 
equilibrium state. In terms of the unphysical rescaled variables we find that the frequencies of the oscillations are 
essentially independent of the value of the coupling constant rj and of the width and amplitude of the Weber- Wheeler 
pulse. We also find that the frequency of X is independent of a while that of P is proportional to yfa. When this 
result is translated back into the physical units we find that the frequency of the scalar field is proportional to the 
mass of the scalar field and the frequency of the vector field is proportional to the mass of the vector field. This result 
is confirmed by investigating two further scenarios. Firstly we consider the evolution of static initial data for the 
string coupled to the gravitational field, but with a Gaussian perturbation to one of the string variables, and secondly 
we consider the same thing but in a Minkowskian background with the gravitational field decoupled. In both cases 
we obtain the same relationship between the frequencies and the mass. 

Having investigated the interaction between a Weber- Wheeler type pulse of gravitational radiation and the cosmic 
string, the next obvious step is to consider the interaction between the string and a pulse of gravitational radiation 
with both polarisation states present. The code will also be a valuable tool in comparing the numerical results with 
those that one gets from using perturbation theory or the thin string limit to approximate the behavior of a cosmic 
string coupled to gravity. 
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